clear 

cd "\\file\Usersw$\wrr15\Home\My Documents\My Files\WU\CHAPTER2"

//log using blogresults, replace

use "replication dta_20211018(ts10 raw).dta"

drop if year>2004 
drop if year<1960

// Create squared variable
gen re_gdpsq = re_gdp^2

// Create missing value variable
gen lemiss = (re_le == .)

rename re_le le
rename re_ts10 ts10
rename re_gdpsq gdpsq
rename re_gdp gdp
rename re_edu edu
rename lre_phealth phealth
rename lre_thealth thealth
rename lre_im im

// Show missing data
mdesc

// Create dummy variables
tab year, gen(year)
tab id, gen (id)

xtset id year
xtreg le ts10 gdp gdpsq edu phealth thealth year2-year45, fe vce(cluster id)

preserve
clear 
use "replication linear imputation.dta"
drop if year>2004 
drop if year<1960
gen re_gdpsq = re_gdp^2
rename re_le le
rename re_ts10 ts10
rename re_gdpsq gdpsq
rename re_gdp gdp
rename re_edu edu
rename lre_phealth phealth
rename lre_thealth thealth
tab year, gen(year)
tab id, gen (id)
xtset id year
xtreg le ts10 gdp gdpsq edu phealth thealth year2-year45, fe vce(cluster id)
restore

// SEM - No FIML
sem (le <- ts10 gdp gdpsq edu phealth thealth id2-id12 year2-year45), vce(cluster id)

/*
// SEM - FIML
// Didn't converge
sem (le <- ts10 gdp gdpsq edu phealth thealth id2-id12 year2-year45), method(mlmv) vce(cluster id)
*/

//Set-up
mi set flong
mi register imputed le ts10 edu phealth thealth im

// Note the default burnin and burnbetween are both 100
// I set "initmcmc(em)" to do 500 iterations because it didn't converge at the default values
mi impute mvn le ts10 edu phealth thealth im  = gdp gdpsq id2-id12 year2-year45, prior(jeffreys) ///
    mcmconly rseed(123) savewlf(wlf, replace)

preserve
use wlf, clear
tsset iter
tsline wlf, ytitle(Worst linear function) xtitle(Burn-in period) name(stable1,replace)
ac wlf, title(Worst linear function) ytitle(Autocorrelations) ///
  note("") name(ac1,replace)
restore

// Stationarity doesn't look great. Increase burnin to 500.

mi impute mvn le ts10 edu phealth thealth im  = gdp gdpsq id2-id12 year2-year45, prior(jeffreys) ///
    mcmconly burnin(500) rseed(123) savewlf(wlf, replace)

preserve
use wlf, clear
tsset iter
tsline wlf, ytitle(Worst linear function) xtitle(Burn-in period) name(stable2,replace)
ac wlf, title(Worst linear function) ytitle(Autocorrelations) ///
  note("") name(ac2,replace)
restore

// Imputation
mi impute mvn le ts10 edu phealth thealth im  = gdp gdpsq id2-id12 year2-year45, prior(jeffreys) ///
    burnin(500) add(10) rseed(123) 

// Analysis
mi xtset id year
mi estimate, cmdok: wcbregress le ts10 gdp gdpsq edu phealth thealth year2-year45, fe vce(cluster id)

how_many_imputations

// Imputation
mi impute mvn le ts10 edu phealth thealth im  = gdp gdpsq id2-id12 year2-year45, prior(jeffreys) ///
    burnin(500) add(`r(add_M)')
	
// Analysis
mi xtset id year
mi estimate: xtreg le ts10 gdp gdpsq edu phealth thealth year2-year45, fe vce(cluster id)

//log close



